library(pacman)
pacman::p_load(tidyverse, brms, ggeffects, kableExtra, tidybayes, install = TRUE)palette = c("#3DC0E1", "#1985A1", "#4C5C68", "#00A896")
palette_group = c("#0C7BDC", "#D55E00")
palette_sd = c(palette[1], palette[3], palette[2])
plot_aes = theme_minimal() +
theme(legend.position = "top",
legend.text = element_text(size = 12),
text = element_text(size = 16, family = "Futura Medium"),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.text = element_text(color = "black"),
axis.line = element_line(colour = "black"),
axis.ticks.y = element_blank())make_table = function(data) {
data %>%
broom.mixed::tidy(conf.int = TRUE) %>%
filter(effect == "fixed") %>%
mutate(term = gsub("\\(Intercept\\)", "intercept", term),
term = gsub("regulation_expression", "signature expression (mindfulness)", term),
term = gsub("active_weekon", "intervention week (on)", term),
term = gsub("active_weekoff", "intervention week (off)", term),
term = gsub("signal_count", "signal", term),
term = gsub("SocialWeekendWeekend", "weekend", term),
term = gsub("gender_f", "gender (female)", term),
term = gsub(":", " x ", term),
`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", estimate, conf.low, conf.high)) %>%
select(term, `b [95% CI]`) %>%
knitr::kable(digits = 2)
}task_data = read.csv("../data/task/cuereact_all062920.csv", stringsAsFactors = FALSE) %>%
filter(grepl("rating", trial_type)) %>%
select(pID, block_num, trial_num, trial_type, resp, rt, stim) %>%
mutate(trial_cond = ifelse(grepl("friend3|friend4", trial_type), "regulation",
ifelse(grepl("friend1|friend2", trial_type), "up-regulation",
ifelse(grepl("mindful", trial_type), "regulation",
ifelse(grepl("nonalc", trial_type), "non-alcohol reactivity",
ifelse(grepl("rating_alc_react", trial_type), "reactivity", NA))))),
trial_cond = factor(trial_cond, levels = c("non-alcohol reactivity", "reactivity", "regulation", "up-regulation"))) %>%
rename("stimulus" = stim) %>%
filter(!is.na(stimulus))
task_conditions = task_data %>%
select(pID, trial_type) %>%
unique() %>%
filter(grepl("friend|mindful", trial_type) | pID %in% c("sub-MURIP012", "sub-MURIP063", "sub-MURIP078")) %>%
mutate(condition = ifelse(grepl("friend", trial_type), "perspective", "mindfulness")) %>%
select(-trial_type) %>%
unique()
task_data_all = task_data %>%
left_join(., task_conditions) %>%
mutate(condition = ifelse(is.na(condition), "control", condition),
trial_cond = factor(trial_cond, c("non-alcohol reactivity", "reactivity", "regulation", "up-regulation")))
task_mindfulness = task_data_all %>%
filter(condition == "mindfulness")file_dir = "../data/dots"
file_pattern = "sub.*"
file_list = list.files(file_dir, pattern = file_pattern)
dots_tmp = data.frame()
for (file in file_list) {
tmp = tryCatch(read.table(file.path(file_dir,file), fill = TRUE, header = TRUE, sep = ",") %>%
extract(file, c("pID", "beta"), ".*(sub-MURIC[0-9]{3}|sub-MURIP[0-9]{3})/(beta_[0-9]{4}).nii") %>%
mutate(stimulus = as.character(stimulus)), error = function(e) message(file))
dots_tmp = rbind(dots_tmp, tmp)
rm(tmp)
}
dots = dots_tmp %>%
left_join(., task_conditions) %>%
mutate(condition = ifelse(is.na(condition), "control", condition)) %>%
group_by(condition) %>%
filter(!stimulus == "nan") %>%
mutate(sd3_signal = 3 * sd(mean_signal, na.rm = TRUE),
grand_mean_signal = mean(mean_signal, na.rm = TRUE),
outlier = ifelse(mean_signal > grand_mean_signal + sd3_signal, 1,
ifelse(mean_signal < grand_mean_signal - sd3_signal, 1, 0))) %>%
select(-sd3_signal, -grand_mean_signal) %>%
filter(!outlier == 1)classifier_data = read.csv("../data/mvpa/logistic_stats.csv", stringsAsFactors = FALSE) %>%
rename("pID" = subjectID) %>%
mutate(predicted_factor = ifelse(predicted > .5, "regulate", "react"),
predicted_factor = as.factor(predicted_factor),
actual_factor = ifelse(actual == 1, "regulate", "react"),
actual_factor = as.factor(actual_factor))merged = task_data_all %>%
full_join(., dots) %>%
filter(!pID == "sub-MURIC303") %>%
filter(trial_cond %in% c("reactivity", "regulation")) %>%
mutate(trial_cond_recode = ifelse(trial_cond == "regulation", .5, -.5),
trial_cond = as.factor(as.character(trial_cond)))
between = merged %>%
select(pID, dot, trial_cond, trial_cond_recode, condition) %>%
group_by(pID, trial_cond, trial_cond_recode, condition) %>%
summarize(dot_between = mean(dot, na.rm = TRUE)) %>%
group_by(trial_cond, condition) %>%
mutate(dot_between_c = scale(dot_between, scale = FALSE, center = TRUE)) %>%
group_by(trial_cond, condition) %>%
mutate(sd_dot = sd(dot_between, na.rm = TRUE),
dot_between_std = dot_between_c / sd_dot,
dot_between_std_noc = dot_between / sd_dot) %>%
select(-sd_dot)
reg_react_expression = between %>%
select(pID, condition, trial_cond, dot_between_std) %>%
unique() %>%
mutate(trial_cond = sprintf("%s_expression", trial_cond)) %>%
spread(trial_cond, dot_between_std)ema = read.csv("../data/ema/cleaned_EMA.csv", stringsAsFactors = FALSE) %>%
select(pID, groupID, Condition, contains("signal"), Session.Name, active_week, drinks_number, contains("Craving"), gender_f, NumberResponses, SocialWeekend, amount_self_c, freq_self_c, contains("React")) %>%
mutate(craving_previous = lag(Craving_Alc),
Condition = gsub("mindful", "mindfulness", Condition),
pID = gsub("muric", "sub-MURIC", pID),
pID = gsub("muri", "sub-MURIP", pID),
active_week = gsub("active", "on", active_week),
active_week = gsub("control", "off", active_week),
time_of_day = ifelse(grepl("Morning", Session.Name), "morning",
ifelse(grepl("Evening", Session.Name), "evening", NA))) %>%
rename("condition" = Condition)
ema_within_intervetion = ema %>%
left_join(., reg_react_expression) %>%
filter(!condition == "control") %>%
group_by(pID) %>%
mutate(drinks_number = scale(drinks_number, center = TRUE, scale = FALSE),
Alc_React_Mindful = scale(Alc_React_Mindful, center = TRUE, scale = FALSE),
Alc_React_Perspective = scale(Alc_React_Perspective, center = TRUE, scale = FALSE),
mindful_scale = scale(Alc_React_Mindful, center = TRUE, scale = TRUE),
perspective_scale = scale(Alc_React_Perspective, center = TRUE, scale = TRUE),
response_scale = ifelse(condition == "mindfulness", mindful_scale, perspective_scale),
craving_previous = scale(craving_previous, center = TRUE, scale = TRUE),
craving = scale(Craving_Alc, center = TRUE, scale = TRUE))
ema_within = ema_within_intervetion %>%
filter(condition == "mindfulness")
ema_between = ema %>%
group_by(pID, condition) %>%
summarize(drinks_between = mean(drinks_number, na.rm = TRUE),
responses_between = mean(Alc_React_Mindful, na.rm = TRUE),
craving_between = mean(Craving_Alc, na.rm = TRUE)) %>%
filter(condition == "mindfulness") %>%
left_join(., reg_react_expression)ema_within %>%
select(pID, drinks_number, mindful_scale, craving, active_week) %>%
gather(key, value, -pID, -active_week) %>%
ggplot(aes(active_week, value, color = key)) +
stat_summary(aes(group = interaction(pID, key)), fun = "mean", geom = "line", color = "grey", size = .2) +
stat_summary(aes(group = key), fun = "mean", geom = "line") +
stat_summary(fun.data = "mean_cl_boot") +
facet_grid(~key) +
plot_aes +
theme(legend.position = "none")Association between how mindful you were when you encountered alcohol and how much you drank (morning signal = reported evening drinking, evening signal = reporting morning/afternoon drinking)
ema_within %>%
ggplot(aes(mindful_scale, drinks_number, color = active_week, fill = active_week)) +
geom_point(size = .5, alpha = .4) +
geom_smooth(aes(group = interaction(pID, active_week)), method = "lm", se = FALSE, size = .2) +
geom_smooth(method = "lm", alpha = .2) +
scale_color_manual(values = palette) +
scale_fill_manual(values = palette) +
plot_aesE.g. association between morning craving and morning/afternoon drinking reported in the evening
ema_within %>%
ggplot(aes(craving_previous, drinks_number, color = active_week, fill = active_week)) +
geom_point(size = .5, alpha = .4) +
geom_smooth(aes(group = interaction(pID, active_week)), method = "lm", se = FALSE, size = .2) +
geom_smooth(method = "lm", alpha = .2) +
scale_color_manual(values = palette) +
scale_fill_manual(values = palette) +
plot_aesE.g. association between morning craving and morning/afternoon mindful responses to alcohol reported in the evening
ema_within %>%
ggplot(aes(craving_previous, mindful_scale, color = active_week, fill = active_week)) +
geom_point(size = .5, alpha = .4) +
geom_smooth(aes(group = interaction(pID, active_week)), method = "lm", se = FALSE, size = .2) +
geom_smooth(method = "lm", alpha = .2) +
scale_color_manual(values = palette) +
scale_fill_manual(values = palette) +
plot_aesE.g. association between morning/afternoon mindful responses to alcohol reported in the evening and craving reported in the evening
ema_within %>%
ggplot(aes(mindful_scale, craving, color = active_week, fill = active_week)) +
geom_point(size = .5, alpha = .4) +
geom_smooth(aes(group = interaction(pID, active_week)), method = "lm", se = FALSE, size = .2) +
geom_smooth(method = "lm", alpha = .2) +
scale_color_manual(values = palette) +
scale_fill_manual(values = palette) +
plot_aesvals = seq(-2, 2, .2)
predicted = ggeffects::ggpredict(drinking_week, c("regulation_expression [vals]", "active_week")) %>%
data.frame()
predicted %>%
ggplot(aes(x, predicted, color = group, fill = group)) +
geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2, color = NA) +
geom_line(size = 1) +
scale_color_manual(name = "intervention week", values = palette[3:4]) +
scale_fill_manual(name = "intervention week", values = palette[3:4]) +
labs(y = "within-person number of drinks\n", x = "\nsignature expression on mindfulness trials (SD)") +
plot_aes +
theme(legend.position = "top")predicted = ggeffects::ggpredict(drinking_week, c("regulation_expression [0, 1]", "active_week")) %>%
data.frame() %>%
mutate(x = ifelse(x == 1, "+1 SD", "mean"),
x = factor(x, levels = c("mean", "+1 SD")),
group = recode(group, "on" = "active", "off" = "control"),
group = factor(group, levels = c("control", "active")))
predicted %>%
filter(x == "mean") %>%
ggplot(aes(group, predicted, color = x, group = x)) +
geom_pointrange(aes(ymin = conf.low, ymax = conf.high), position = position_dodge(0.05), size = 1.5) +
geom_line(position = position_dodge(0.05), size = 1.5) +
scale_color_manual(name = "mindfulness signature expression", values = palette) +
scale_x_discrete(expand = c(.1, .1)) +
coord_cartesian(ylim = c(-.3, .3)) +
labs(y = "within-person number of drinks\n", x = "\nintervention week") +
plot_aes +
theme(legend.position = "top")predicted %>%
ggplot(aes(group, predicted, color = x, group = x)) +
geom_pointrange(aes(ymin = conf.low, ymax = conf.high), position = position_dodge(0.05), size = 1.5) +
geom_line(position = position_dodge(0.05), size = 1.5) +
scale_color_manual(name = "mindfulness signature expression", values = palette) +
scale_x_discrete(expand = c(.1, .1)) +
coord_cartesian(ylim = c(-.3, .3)) +
labs(y = "within-person number of drinks\n", x = "\nintervention week") +
plot_aes +
theme(legend.position = "top")make_table(drinking_week)| term | b [95% CI] |
|---|---|
| intercept | 0.20 [0.06, 0.34] |
| intervention week (on) | -0.07 [-0.22, 0.08] |
| signature expression (mindfulness) | 0.11 [0.01, 0.21] |
| signal | -0.01 [-0.01, -0.00] |
| intervention week (on) x signature expression (mindfulness) | -0.23 [-0.38, -0.08] |
hypothesis(drinking_week,
c(expression_mean = "active_weekon = 0",
expression_1sd = "active_weekon + active_weekon:regulation_expression = 0"))$hypothesis %>%
mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", Estimate, CI.Lower, CI.Upper)) %>%
select(Hypothesis, `b [95% CI]`) %>%
knitr::kable()| Hypothesis | b [95% CI] |
|---|---|
| expression_mean | -0.07 [-0.22, 0.08] |
| expression_1sd | -0.30 [-0.50, -0.09] |
vals = seq(-2, 2, .2)
predicted = ggeffects::ggpredict(response_week, c("regulation_expression [vals]", "active_week")) %>%
data.frame()
predicted %>%
ggplot(aes(x, predicted, color = group, fill = group)) +
geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .2, color = NA) +
geom_line(size = 1) +
scale_color_manual(name = "intervention week", values = palette[3:4]) +
scale_fill_manual(name = "intervention week", values = palette[3:4]) +
labs(y = "within-person mindful response rating (SD)\n", x = "\nsignature expression on mindfulness trials (SD)") +
plot_aes +
theme(legend.position = "top")predicted = ggeffects::ggpredict(response_week, c("regulation_expression [0, 1]", "active_week")) %>%
data.frame() %>%
mutate(x = ifelse(x == 1, "+1 SD", "mean"),
x = factor(x, levels = c("mean", "+1 SD")),
group = recode(group, "on" = "active", "off" = "control"),
group = factor(group, levels = c("control", "active")))
predicted %>%
filter(x == "mean") %>%
ggplot(aes(group, predicted, color = x, group = x)) +
geom_pointrange(aes(ymin = conf.low, ymax = conf.high), position = position_dodge(0.05), size = 1.5) +
geom_line(position = position_dodge(0.05), size = 1.5) +
scale_color_manual(name = "mindfulness signature expression", values = palette) +
scale_x_discrete(expand = c(.1, .1)) +
coord_cartesian(ylim = c(-.5, .65)) +
labs(y = "within-person mindful response rating (SD)\n", x = "\nintervention week") +
plot_aes +
theme(legend.position = "top")predicted %>%
ggplot(aes(group, predicted, color = x, group = x)) +
geom_pointrange(aes(ymin = conf.low, ymax = conf.high), position = position_dodge(0.05), size = 1.5) +
geom_line(position = position_dodge(0.05), size = 1.5) +
scale_color_manual(name = "mindfulness signature expression", values = palette) +
scale_x_discrete(expand = c(.1, .1)) +
coord_cartesian(ylim = c(-.5, .65)) +
labs(y = "within-person mindful response rating (SD)\n", x = "\nintervention week") +
plot_aes +
theme(legend.position = "top")make_table(response_week)| term | b [95% CI] |
|---|---|
| intercept | -0.49 [-0.73, -0.27] |
| intervention week (on) | 0.44 [0.19, 0.69] |
| signature expression (mindfulness) | -0.10 [-0.24, 0.06] |
| signal | 0.01 [0.00, 0.02] |
| intervention week (on) x signature expression (mindfulness) | 0.26 [-0.01, 0.51] |
hypothesis(response_week,
c(expression_mean = "active_weekon = 0",
expression_1sd = "active_weekon + active_weekon:regulation_expression = 0"))$hypothesis %>%
mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", Estimate, CI.Lower, CI.Upper)) %>%
select(Hypothesis, `b [95% CI]`) %>%
knitr::kable()| Hypothesis | b [95% CI] |
|---|---|
| expression_mean | 0.44 [0.19, 0.69] |
| expression_1sd | 0.71 [0.36, 1.03] |
hypothesis(
fit_brm,
'b_drinksnumber_mindful_scale * b_mindfulscale_active_weekon + cor_pID__mindfulscale_active_weekon__drinksnumber_active_weekon * sd_pID__mindfulscale_active_weekon * sd_pID__drinksnumber_mindful_scale = 0',
class = NULL,
seed = 6523
)## Hypothesis Tests for class :
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (b_drinksnumber_m... = 0 -0.31 0.14 -0.59 -0.03 NA
## Post.Prob Star
## 1 NA *
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.